Cross-tissue by OTD
check correlations
ggplot(data,aes(x=pve50,y=pge50,ymin=pge025,ymax=pge975,col=`pge025>0.01`)) + geom_pointrange(col='gray') + geom_point()

cor.test(hun$pge50,hun$pve50)
##
## Pearson's product-moment correlation
##
## data: hun$pge50 and hun$pve50
## t = 130.5, df = 16948, p-value < 2.2e-16
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## 0.7003744 0.7153942
## sample estimates:
## cor
## 0.7079643
ggplot(data,aes(x=log10(n_gamma50),y=pge50,col=`pge025>0.01`)) + geom_point(alpha=0.3)

check tissue-specific correlations PGE vs PVE
tislist <- scan('/Volumes/im-lab/nas40t2/hwheeler/cross-tissue/nine.tissue.list',sep="\n",what="character")
ts <- data.frame()
for(tis in tislist){
data <- read.table(my.vol %&% tis %&% '_TS_exp_BSLMM-s100K_iterations_all_chr1-22_2015-08-06.txt',header=T,sep="\t")
subdata <- select(data,gene,pve50,pge50,pge025,pge975) %>% mutate(tissue=tis,`pge025>0.01`=pge025>0.01)
res<-cor.test(subdata$pge50,subdata$pve50)
cat(tis,"\tPearson R=",round(res$estimate,3),"\tP-value=",res$p.value,"\n")
ts <- rbind(ts,subdata)
}
## Adipose-Subcutaneous Pearson R= 0.441 P-value= 0
## Artery-Tibial Pearson R= 0.473 P-value= 0
## Heart-LeftVentricle Pearson R= 0.443 P-value= 0
## Lung Pearson R= 0.426 P-value= 0
## Muscle-Skeletal Pearson R= 0.537 P-value= 0
## Nerve-Tibial Pearson R= 0.518 P-value= 0
## Skin-SunExposed(Lowerleg) Pearson R= 0.501 P-value= 0
## Thyroid Pearson R= 0.572 P-value= 0
## WholeBlood Pearson R= 0.558 P-value= 0
##exclude gene with pve ~0 and pge =1
#ts<-filter(ts,gene!="ENSG00000196409.7")
ggplot(ts,aes(x=pve50,y=pge50,ymin=pge025,ymax=pge975,color=`pge025>0.01`) ) + facet_wrap(~tissue,ncol=3) + geom_pointrange(col='gray')+geom_point()+theme_bw()+ggtitle("BSLMM Tissue-Specific")

check tissue-wide correlations PGE v PVE
tislist <- scan('/Volumes/im-lab/nas40t2/hwheeler/cross-tissue/nine.tissue.list',sep="\n",what="character")
tw <- data.frame()
for(tis in tislist){
data <- read.table(my.vol %&% tis %&% '_TW_exp_BSLMM-s100K_iterations_all_chr1-22_2015-10-18.txt',header=T,sep="\t")
subdata <- select(data,gene,pve50,pge50,pge025,pge975) %>% mutate(tissue=tis,`pge025>0.01`=pge025>0.01)
res<-cor.test(subdata$pge50,subdata$pve50)
cat(tis,"\tPearson R=",round(res$estimate,3),"\tP-value=",res$p.value,"\n")
tw <- rbind(tw,subdata)
}
## Adipose-Subcutaneous Pearson R= 0.653 P-value= 0
## Artery-Tibial Pearson R= 0.66 P-value= 0
## Heart-LeftVentricle Pearson R= 0.589 P-value= 0
## Lung Pearson R= 0.611 P-value= 0
## Muscle-Skeletal Pearson R= 0.662 P-value= 0
## Nerve-Tibial Pearson R= 0.671 P-value= 0
## Skin-SunExposed(Lowerleg) Pearson R= 0.654 P-value= 0
## Thyroid Pearson R= 0.67 P-value= 0
## WholeBlood Pearson R= 0.663 P-value= 0
ggplot(tw,aes(x=pve50,y=pge50,ymin=pge025,ymax=pge975,color=`pge025>0.01`) ) + facet_wrap(~tissue,ncol=3) + geom_pointrange(col='gray')+geom_point()+theme_bw()+ggtitle("BSLMM Tissue-Wide")

plot tissue-specific PVE vs GCTA marginal h2
tislist <- scan('/Volumes/im-lab/nas40t2/hwheeler/cross-tissue/nine.tissue.list',sep="\n",what="character")
tislist2 <- scan('/Volumes/im-lab/nas40t2/hwheeler/cross-tissue/nine.spaces.tissue.list',sep="\n",what="character")
ts <- data.frame()
for(i in 1:length(tislist)){
bs <- read.table(my.vol %&% tislist[i] %&% '_TS_exp_BSLMM-s100K_iterations_all_chr1-22_2015-08-06.txt',header=T,sep="\t") %>% select(gene,pve50)
h2 <- read.table("/Volumes/im-lab/nas40t2/hwheeler/cross-tissue/gtex-h2-estimates/GTEx.resid.tissue-specific.h2_" %&% tislist2[i] %&% "_marginal.local_2015-03-23.txt",header=T, sep="\t") %>% select(tissue,ensid,h2) %>% mutate(gene=ensid)
subdata <- inner_join(h2,bs,by="gene")
res<-cor.test(subdata$pve50,subdata$h2)
cat(tislist2[i],"\tPearson R=",round(res$estimate,3),"\tP-value=",res$p.value,"\n")
ts <- rbind(ts,subdata)
}
## Adipose - Subcutaneous Pearson R= 0.565 P-value= 0
## Artery - Tibial Pearson R= 0.611 P-value= 0
## Heart - Left Ventricle Pearson R= 0.589 P-value= 0
## Lung Pearson R= 0.578 P-value= 0
## Muscle - Skeletal Pearson R= 0.629 P-value= 0
## Nerve - Tibial Pearson R= 0.587 P-value= 0
## Skin - Sun Exposed (Lower leg) Pearson R= 0.648 P-value= 0
## Thyroid Pearson R= 0.608 P-value= 0
## Whole Blood Pearson R= 0.688 P-value= 0
ggplot(ts,aes(x=h2,y=pve50))+geom_point(alpha=0.4)+coord_cartesian(xlim=c(0,1),ylim=c(0,1))+xlab(expression("GCTA h"^2))+ylab('BSLMM PVE')+geom_abline(c(0,1),color='red') + facet_wrap(~tissue,ncol=3)+theme_bw()+ggtitle("Tissue-Specific")

plot tissue-wide PVE vs GCTA marginal h2
tislist <- scan('/Volumes/im-lab/nas40t2/hwheeler/cross-tissue/nine.tissue.list',sep="\n",what="character")
tw <- data.frame()
for(i in 1:length(tislist)){
bs <- read.table(my.vol %&% tislist[i] %&% '_TW_exp_BSLMM-s100K_iterations_all_chr1-22_2015-10-18.txt',header=T,sep="\t") %>% select(gene,pve50)
h2 <- read.table("/Volumes/im-lab/nas40t2/hwheeler/cross-tissue/gtex-h2-estimates/GTEx.tissue-wide.h2_" %&% tislist[i] %&% "_marginal.local_2015-03-24.txt",header=T, sep="\t") %>% select(tissue,ensid,h2) %>% mutate(gene=ensid)
subdata <- inner_join(h2,bs,by="gene")
res<-cor.test(subdata$pve50,subdata$h2)
cat(tislist[i],"\tPearson R=",round(res$estimate,3),"\tP-value=",res$p.value,"\n")
tw <- rbind(tw,subdata)
}
## Adipose-Subcutaneous Pearson R= 0.414 P-value= 0
## Artery-Tibial Pearson R= 0.414 P-value= 0
## Heart-LeftVentricle Pearson R= 0.307 P-value= 0
## Lung Pearson R= 0.377 P-value= 0
## Muscle-Skeletal Pearson R= 0.459 P-value= 0
## Nerve-Tibial Pearson R= 0.442 P-value= 0
## Skin-SunExposed(Lowerleg) Pearson R= 0.459 P-value= 0
## Thyroid Pearson R= 0.465 P-value= 0
## WholeBlood Pearson R= 0.38 P-value= 0
ggplot(tw,aes(x=h2,y=pve50))+geom_point(alpha=0.4)+coord_cartesian(xlim=c(0,1),ylim=c(0,1))+xlab(expression("GCTA h"^2))+ylab('BSLMM PVE')+geom_abline(c(0,1),color='red') + facet_wrap(~tissue,ncol=3)+theme_bw()+ggtitle("Tissue-Wide")
